library(magrittr)
library(tidyverse)
library(Seurat)
library(readxl)
library(cowplot)
library(colorblindr)
library(viridis)
library(progeny)
library(destiny)

coi <- params$cell_type_super
cell_sort <- params$cell_sort
cell_type_major <- params$cell_type_major
louvain_resolution <- params$louvain_resolution
louvain_cluster <- params$louvain_cluster

1 Cluster markers

1.1 Major T.super markers for cell assign

### load all data ---------------------------------
source("_src/global_vars.R")

# seu_obj <- read_rds(paste0("/work/shah/isabl_data_lake/analyses/16/52/1652/celltypes/", coi, "_processed.rds"))
seu_obj <- read_rds(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/outs_pre/", coi, "_seurat_", louvain_resolution, ".rds"))
# seu_obj <- read_rds(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_highqc.rds"))

myfeatures <- c("umapharmony_1", "umapharmony_2", "sample", louvain_cluster, "doublet", "nCount_RNA", "nFeature_RNA", "percent.mt", "doublet_score")

plot_data_wrapper <- function(cluster_res) {
  cluster_res <- enquo(cluster_res)
  as_tibble(FetchData(seu_obj, myfeatures)) %>% 
    left_join(meta_tbl, by = "sample") %>% 
    rename(cluster = !!cluster_res) %>% 
    mutate(cluster = as.character(cluster),
           tumor_supersite = ordered(tumor_supersite, levels = names(clrs$tumor_supersite)))
}

plot_data <- plot_data_wrapper(louvain_cluster)

1.2 Subtype currated markers

helper_f <- function(x) ifelse(is.na(x), "", x)

markers_v6_super[[coi]] %>% 
  group_by(subtype) %>% 
  mutate(rank = row_number(gene)) %>% 
  spread(subtype, gene) %>% 
  mutate_all(.funs = helper_f) %>% 
  formattable::formattable()
rank CD4.T.CXCL13 CD4.T.naive CD4.T.reg CD4.Th17 CD8.T.CXCL13 CD8.T.cytotoxic CD8.T.ISG Cycling.T.NK dissociated NK.CD56 NK.cytotoxic
1 CD4 CCR7 CD4 CCR6 CCL3 CCL4 IFI6 ASPM BTG1 AREG FCGR3A
2 CD40LG CD4 FOXP3 IL4I1 CCL4L2 CCL5 IFIT1 CENPF DNAJB1 FCER1G FGFBP2
3 CTLA4 CD40LG IL2RA IL7R CD8A CD8A IFIT2 HIST1H4C DUSP1 GNLY KLRD1
4 CXCL13 IL7R TNFRSF4 KLRB1 CRTAM CD8B IFIT3 HMGB2 EGR1 KLRC1 KLRF1
5 FKBP5 KLF2 TRAC LST1 CXCL13 CRTAM ISG15 MKI67 FOS KRT81 PRF1
6 IL6ST LTB LTB FABP5 CST7 MX1 STMN1 FOSB TRDC SPON2
7 ITM2A TCF7 TNFSF13B GZMB DTHD1 MX2 TOP2A HSPA1A TYROBP
8 MAF TPT1 HAVCR2 GZMA RSAD2 TUBA1B HSPA1B XCL1
9 NMB IFNG GZMH TUBB HSPA6 XCL2
10 NR3C1 LAG3 GZMK TYMS JUN
11 PDCD1 MIR155HG GZMM JUNB
12 TNFRSF4 PHLDA1 HLA-DPB1 KLF2
13 TOX2 PTMS ITM2C MT1E
14 TSHZ2 RBPJ KLRG1 MT1X
15 TNFRSF9 TRGC2

1.3 Subtype cluster markers

## define patient specific clusters
patient_specific_clusters <- plot_data %>% 
  group_by(patient_id_short, cluster) %>% 
  tally %>% 
  group_by(cluster) %>% 
  mutate(nrel = n / sum(n),
         ntotal = sum(n)) %>% 
  arrange(desc(nrel)) %>% 
  distinct(cluster, .keep_all = T) %>% 
  filter(nrel > 0.5) %>% 
  ungroup %>% 
  mutate(cluster_label_ps = make.names(paste0("OV.", patient_id_short, ".specific"), unique = T),
         cluster = as.numeric(cluster)) %>% 
  distinct(cluster, cluster_label_ps)



## Hypergeometric test --------------------------------------

# marker_tbl <- read_tsv(paste0("/work/shah/isabl_data_lake/analyses/16/52/1652/celltypes/", coi, "_markers.tsv")) %>% 
#   filter(resolution == louvain_resolution)
marker_tbl <- read_tsv(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/outs_pre/", coi, "_markers_", louvain_resolution, ".tsv"))
# marker_tbl <- read_tsv(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_highqc_markers_02.tsv"))

test_set <- marker_tbl %>% 
  group_by(cluster) %>% 
  filter(!str_detect(gene, "^RPS|^RPL")) %>% 
  slice(1:50) %>% 
  mutate(k = length(cluster)) %>% 
  ungroup %>%
  select(cluster, gene, k) %>% 
  mutate(join_helper = 1) %>% 
  group_by(cluster, join_helper, k) %>%
  nest(test_set = gene)

markers_doub_tbl <- markers_v6 %>% 
  enframe("subtype", "gene") %>% 
  filter(!(subtype %in% unique(c(coi, cell_type_major)))) %>% 
  unnest(gene) %>% 
  group_by(gene) %>% 
  filter(length(gene) == 1) %>% 
  mutate(subtype = paste0("doublet.", subtype)) %>% 
  bind_rows(tibble(subtype = "Mito.high", gene = grep("^MT-", rownames(seu_obj), value = T)))

ref_set <- markers_v6_super[[coi]] %>% 
  bind_rows(markers_doub_tbl) %>% 
  group_by(subtype) %>% 
  mutate(m = length(gene),
         n = length(rownames(seu_obj))-m,
         join_helper = 1) %>% 
  group_by(subtype, m, n, join_helper) %>%
  nest(ref_set = gene)

hyper_tbl <- test_set %>% 
  left_join(ref_set, by = "join_helper") %>% 
  group_by(cluster, subtype, m, n, k) %>%
  do(q = length(intersect(unlist(.$ref_set), unlist(.$test_set)))) %>%
  mutate(pval = 1-phyper(q = q, m = m, n = n, k = k)) %>%
  ungroup %>%
  mutate(qval = p.adjust(pval, "BH"),
         sig = qval < 0.01)

# hyper_tbl %>% 
#   group_by(subtype) %>% 
#   filter(any(qval < 0.01)) %>%
#   ggplot(aes(subtype, -log10(qval), fill = sig)) +
#   geom_bar(stat = "identity") +
#   facet_wrap(~cluster) +
#   coord_flip()
  
low_rank <- str_detect(unique(hyper_tbl$subtype), "doublet|dissociated")
subtype_lvl <- c(sort(unique(hyper_tbl$subtype)[!low_rank]), sort(unique(hyper_tbl$subtype)[low_rank]))
  
cluster_label_tbl <- hyper_tbl %>% 
  mutate(subtype = ordered(subtype, levels = subtype_lvl)) %>% 
  arrange(qval, subtype) %>%
  group_by(cluster) %>% 
  slice(1) %>% 
  mutate(subtype = ifelse(sig, as.character(subtype), paste0("unknown_", cluster))) %>% 
  select(cluster, cluster_label = subtype) %>% 
  ungroup %>% 
  mutate(cluster_label = make.unique(cluster_label, sep = "_")) %>% 
  left_join(patient_specific_clusters, by = "cluster") %>% 
  mutate(cluster_label = ifelse(is.na(cluster_label_ps), cluster_label, cluster_label_ps)) %>% 
  select(-cluster_label_ps)

seu_obj$cluster_label <- unname(deframe(cluster_label_tbl)[as.character(unlist(seu_obj[[paste0("RNA_snn_res.", louvain_resolution)]]))])
plot_data$cluster_label <- seu_obj$cluster_label

cluster_n_tbl <- seu_obj$cluster_label %>% 
  table() %>% 
  enframe("cluster_label", "cluster_n") %>% 
  mutate(cluster_nrel = cluster_n/sum(cluster_n))

marker_sheet <- marker_tbl %>% 
  left_join(cluster_label_tbl, by = "cluster") %>% 
  group_by(cluster_label) %>% 
  filter(!str_detect(gene, "^RPS|^RPL")) %>% 
  slice(1:50) %>% 
  mutate(rank = row_number(-avg_logFC)) %>% 
  select(cluster_label, gene, rank) %>% 
  ungroup %>% 
  mutate(cluster_label = ordered(cluster_label, levels = unique(c(names(clrs$cluster_label[[coi]]), sort(cluster_label))))) %>% 
  spread(cluster_label, gene) %>% 
  mutate_all(.funs = helper_f)

marker_tbl_annotated <- marker_tbl %>% 
  left_join(cluster_label_tbl, by = "cluster") %>% 
  left_join(cluster_n_tbl, by = "cluster_label") %>% 
  select(-cluster, -resolution) %>% 
  mutate(cluster_label = ordered(cluster_label, levels = unique(c(names(clrs$cluster_label[[coi]]), sort(cluster_label))))) %>% 
  arrange(cluster_label, -avg_logFC, p_val_adj)

write_tsv(marker_sheet, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_marker_sheet.tsv"))

write_tsv(marker_sheet, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/supplementary_tables/", coi, "_marker_sheet.tsv"))

write_tsv(marker_tbl_annotated, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_marker_table_annotated.tsv"))

write_tsv(marker_tbl_annotated, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/supplementary_tables/", coi, "_marker_table_annotated.tsv"))

formattable::formattable(marker_sheet)
rank CD4.T.naive CD4.T.CXCL13 CD4.T.reg CD4.Th17 CD8.T.cytotoxic CD8.T.CXCL13 CD8.T.ISG NK.CD56 NK.cytotoxic Cycling.T.NK Cycling.T.NK_1 dissociated dissociated_1 doublet.Fibroblast doublet.Monocyte OV.008.specific OV.008.specific.1 OV.085.specific
1 IL7R CXCL13 TNFRSF4 KLRB1 GZMK CXCL13 IFIT3 GNLY FGFBP2 STMN1 CENPF CCL4L2 HSPA1A DCN HLA-DRA KLF2 SOX4 IGFBP5
2 CCR7 NMB IL2RA IL4I1 CD8A GZMB ISG15 TYROBP FCGR3A TUBA1B ASPM CCL4 HSPA1B IGFBP5 CST3 CCR7 PTCRA TAGLN
3 KLF2 NR3C1 FOXP3 IL7R CD8B CCL4L2 MX1 AREG SPON2 MKI67 MKI67 IFNG MT1X RBP1 CXCL8 JUNB MAL ADIRF
4 EEF1B2 MAF CTLA4 LTB ITM2C MIR155HG IFIT1 KLRC1 PRF1 HIST1H4C HMGB2 FOS DNAJB1 C7 SPP1 FOS MZB1 DCN
5 TPT1 FKBP5 LTB LST1 GZMH TNFRSF9 RSAD2 FCER1G KLRF1 TUBB TOP2A FOSB MT1E MEG3 S100A9 DUSP1 DNTT IGFBP7
6 EEF1A1 IL6ST RTKN2 TNFSF13B CCL5 HAVCR2 IFIT2 TRDC GNLY TOP2A UBE2C TNF HSPA6 CALD1 S100A8 SELL TFDP2 SPARCL1
7 MAL ITM2A BATF CCR6 TRGC2 RBPJ IFI6 XCL1 KLRD1 TYMS CCNB1 JUN HSP90AA1 IGFBP4 HLA-DQA1 IL7R CD1E CALD1
8 TCF7 TSHZ2 TNFRSF18 CTSH KLRG1 LAG3 MX2 KLRD1 NKG7 CENPF UBE2S EGR1 HSPH1 RARRES2 CD74 AREG STMN1 MGP
9 CD40LG CTLA4 SAT1 AQP3 GZMA IFNG IFI44L KRT81 CX3CR1 HMGB2 PTTG1 CCL3L1 HSPE1 NR2F2 LYZ EEF1A1 ARPP21 C11orf96
10 SELL CD40LG TBC1D4 CCL20 CCL4 CCL3 ISG20 XCL2 GZMB ASPM TPX2 NFKBID HSPD1 SELENOP FTL CD69 AC084033.3 MYL9
11 GPR183 PDCD1 TIGIT NFKBIA CRTAM PTMS HERC5 IGFBP2 PLAC8 NUSAP1 STMN1 CD69 HSPB1 MDK APOE BTG2 CDK6 TIMP3
12 LDHB CD4 GADD45A RORA CST7 CD8A SAMD9L CLIC3 CLIC3 PCLAF KPNA2 NR4A2 HSP90AB1 SOX4 HLA-DQB1 GPR183 MAP1A MDK
13 NOSIP LIMS1 TNFRSF1B TNFRSF25 GZMM CRTAM OAS1 IL2RB PLEK HMGN2 CENPE AC020916.1 DNAJA1 ADIRF MARCKS PIK3IP1 AC011893.1 TPM2
14 SNHG8 TNFRSF4 PMAIP1 CEBPD HLA-DPB1 PHLDA1 TNFSF10 CEBPD TYROBP H2AFZ CKS2 EGR2 JUN MGP AIF1 CD55 GLUL NR2F2
15 PABPC1 RNF19A UGP2 TNFAIP3 DTHD1 FABP5 STAT1 KRT86 PTGDS PCNA TUBA1B DUSP1 CACYBP EGR1 C1QB DNAJB1 ADA FBLN1
16 NOP53 CORO1B IKZF2 NCR3 PPP1R14B TIGIT EIF2AK2 TXK EFHD2 HIST1H1B HMMR TNFSF9 HSPA8 STAR FTH1 FKBP5 AL357060.1 IGF1
17 LEF1 RBPJ TNFRSF9 SLC4A10 CD3G KRT86 OAS3 CTSW FCER1G DUT DLGAP5 KLF6 FKBP4 SFRP4 BASP1 TSC22D3 GRASP RAMP1
18 EIF3E CPM ICOS TMIGD2 THEMIS JAML SAMD9 KLRB1 CST7 CLSPN CCNB2 TNFAIP3 CHORDC1 CLU C1QA RACK1 CD1B SELENOP
19 LTB ZBED2 LINC01943 DPP4 DUSP2 CCL5 XAF1 MATK GZMH SMC4 TUBB4B DUSP2 RGS2 ADAMTS1 C15orf48 LDHB MIR181A1HG RARRES2
20 RACK1 AC004585.1 IL32 TPT1 CD3D LINC01871 GBP1 CCL3 ADGRG1 ATAD2 NUSAP1 BTG2 FOS TIMP2 FN1 SARAF CCDC26 LUM
21 NACA TOX2 SOX4 MYBL1 LYAR CXCR6 EPSTI1 CD7 CCL3 SMC2 BIRC5 IER2 ANXA1 TCEAL4 MNDA CXCR4 JCHAIN IGFBP6
22 UBA52 DUSP4 ARID5B S100A4 TC2N TNIP3 IFI44 NKG7 HOPX TPX2 HMGN2 PPP1R15A DNAJB4 C1R G0S2 PLAC8 VIPR2 CARMN
23 TOMM7 AHI1 CD27 CD40LG EOMES PDCD1 PLSCR1 CD63 IGFBP7 TMPO ARL6IP1 NR4A1 PPP1R15A WFDC2 APOC1 EEF1B2 ID1 IGFBP4
24 SOCS3 ICA1 BIRC3 SPOCK2 CXCR6 HLA-DRB1 IFI35 HOPX ZEB2 HELLS CDC20 NFKBIA ZFAND2A FHL2 SOD2 TCF7 SOCS2 DST
25 JUNB ARID5B LAYN LINC01871 HLA-DPA1 GAPDH USP18 TMIGD2 PRSS23 UBE2C CDKN3 JUND FOSB IFITM3 NPC2 TPT1 RCAN1 CAV1
26 SERINC5 CD84 CORO1B ERN1 HLA-DRB1 FAM3C OASL CMC1 AKR1C3 NASP SMC4 GADD45B DNAJA4 PEG3 LST1 FOSB CLDN5 ID4
27 TMEM123 CCDC50 TYMP JAML SLF1 CTLA4 MT2A TNFRSF18 CD247 DEK CKS1B TSC22D3 DUSP1 C1S GSN NOSIP CASC15 NUPR1
28 EEF2 IGFL2 DUSP4 FKBP11 APOBEC3G SPRY1 HELZ2 KLRC2 MYBL1 RRM2 KIF20B TAGAP TSPYL2 LUM MS4A6A SC5D PFKFB2 CSRP2
29 FXYD5 RGS1 CD4 IFNGR1 PPP2R5C CCND2 TRIM22 SRGAP3 AREG CKS1B GTSE1 ZFP36 SERPINH1 SERPINF1 IL1B NACA GALNT2 CDKN1C
30 TSHZ2 BATF ENTPD1 MGAT4A KIAA1551 CD63 CMPK2 GSTP1 C1orf21 KNL1 TUBA1C ZFP36L1 UBC CEBPD PSAP AP3M2 HES4 PGR
31 TRABD2A SRGN CTSC S100A6 F2R CD8B LY6E LAT2 S1PR5 HMGB1 TUBB RGCC UBB AKAP12 GRN BTG1 MARCKSL1 COL6A1
32 ANK3 CH25H MIR4435-2HG ELK3 SLAMF7 GZMH PARP14 GZMB CTSW MCM7 SGO2 IL7R AHSA1 TIMP1 CD83 NOP53 TP53INP1 COL6A2
33 SARAF SPOCK2 LINC02099 EEF1A1 CXCR4 ENTPD1 NT5C3A LINC00996 ABHD17A FABP5 H2AFZ CRTAM NEU1 CST3 CTSH ZBTB16 APBA2 EMX2
34 AQP3 ZNRF1 MAGEH1 KIT SH2D1A SRGAP3 DDX58 PRF1 KLF2 UBE2S JPT1 KDM6B DNAJB6 NR2F1 GLUL ERAP2 TSHR PALLD
35 RIPOR2 CHN1 SPOCK2 LTC4S FAM102A VCAM1 IRF7 KLRF1 PTPN12 GAPDH CEP55 NR4A3 GADD45B DLK1 MEF2C CCND3 NREP PPP1R14A
36 AP3M2 TNFRSF25 PHACTR2 RUNX2 CD52 ID2 RNF213 NCAM1 TTC38 CXCL13 NUF2 ATF3 KLF6 SERPING1 CYBB EEF1D MME RBP1
37 ZFAS1 CD200 CARD16 RORC YBX3 HLA-DRA DDX60L CXXC5 KLRB1 RANBP1 KIF14 DDX3X BTG2 BEX3 CTSB PPP1R15A AC002454.1 C7
38 LINC02273 METTL8 S100A4 ZBTB16 STK17A ITGAE DDX60 SH2D1B CCL4 H2AFX KNL1 DUSP6 JUNB C11orf96 CSF3R PLK3 CHI3L2 MFGE8
39 EIF4B RILPL2 STAM FAM241A CCR5 DUSP4 SAT1 MCTP2 PTGDR MCM3 CENPA GZMK TNF FILIP1L CD14 ZFP36L2 SMIM3 KANK2
40 TOB1 TNFRSF18 GLRX PDE4D GPR174 LYST PPM1K IFITM3 ITGB2 EZH2 CDK1 KLRG1 ERN1 RARRES1 C1QC EEF2 SSBP2 NR2F1
41 SESN3 SLA SPATS2L IL23R COTL1 TNFSF4 PNPT1 ZNF683 XBP1 TUBB4B HMGB3 JUNB JUND COL1A2 SGK1 HNRNPA1 UHRF1 MFAP4
42 NSA2 SMCO4 AC005224.3 B3GALT2 CCL4L2 NDFIP2 PARP9 ITGA1 CEP78 H2AFV CDCA8 RASGEF1B ATF3 NUPR1 SPI1 VSIR LRRC28 COL1A1
43 PASK BTLA MAF CERK CD3E AKAP5 IFIH1 IFITM2 ARL4C CDK1 PLK1 ANXA1 DEDD2 TCEAL9 FCGRT PASK BCL11A PBX1
44 TNFRSF25 NAP1L4 PBXIP1 TLE1 TUBA4A GOLIM4 SP110 CCL5 BIN2 HIST1H1D AURKA MCL1 CD69 MARCKSL1 EGR1 EIF3H SCAI GSN
45 FAU FYB1 F5 EEF1B2 PECAM1 CD27 OAS2 ITGAX LITAF HNRNPAB TROAP CXCR4 CLK1 SLC40A1 FCGR2A TXNIP ATP6AP1L PDGFRB
46 EEF1D MIR155HG SLAMF1 CFH ARAP2 SNAP47 C19orf66 CD38 TRDC CENPE H2AFV IER5 IER5L CFH PLAUR LDLRAP1 RUFY3 SERPINF1
47 LDLRAP1 PTPN13 BTG3 PERP ITGA1 RGS1 STAT2 SAMD3 TXK TK1 KIF2C MYADM H3F3B APOE CPVL CMTM8 CD79A SFRP1
48 CTSL SESN3 TRAC PLAT JAML HLA-DPA1 LAG3 SLC16A3 MYOM2 ZWINT MAD2L1 CD8A NR4A1 GNG11 MS4A7 SCML1 GNA15 LHFPL6
49 ITGA6 BIRC3 IL1R1 KIF5C CD84 LINC02446 LAP3 CAPN12 GZMM BIRC5 NUCKS1 PTGER4 CXCR4 CDKN1C ALDH2 LINC00402 HHIP-AS1 PLAC9
50 PFDN5 TP53INP1 DNPH1 PLCB1 AOAH SAMSN1 APOL6 CD247 CD300A PTTG1 RAD21 ZFP36L2 EIF4A2 NBL1 SERPINA1 RIPK2 GSTM3 SERPING1

2 Clusters

2.1 sizes

enframe(sort(table(seu_obj$cluster_label))) %>% 
  mutate(name = ordered(name, levels = rev(name))) %>% 
  ggplot() +
  geom_bar(aes(name, value), stat = "identity") +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5)) +
  labs(y = c("#cells"), x = "cluster")

2.2 UMAP

alpha_lvl <- ifelse(nrow(plot_data) < 20000, 0.2, 0.1)
pt_size <- ifelse(nrow(plot_data) < 20000, 0.2, 0.05)

common_layers_disc <- list(  
  geom_point(size = pt_size, alpha = alpha_lvl),
  NoAxes(),
  guides(color = guide_legend(override.aes = list(size = 2, alpha = 1))),
  labs(color = "")
)

common_layers_cont <- list(  
  geom_point(size = pt_size, alpha = alpha_lvl),
  NoAxes(),
  scale_color_gradientn(colors = viridis(9)),
  guides(color = guide_colorbar())
)

ggplot(plot_data, aes(umapharmony_1, umapharmony_2, color = cluster_label)) + 
  common_layers_disc +
  #facet_wrap(~therapy) +
  ggtitle("Sub cluster")

3 Filtering out doublet clusters

3.1 UMAP

my_subtypes <- names(clrs$cluster_label[[coi]])
my_subtypes <- c(my_subtypes, unlist(lapply(paste0("_", 1:3), function(x) paste0(my_subtypes, x)))) %>% .[!str_detect(., "doublet|dissociated")]
my_subtypes <- my_subtypes[my_subtypes %in% unique(seu_obj$cluster_label)]
my_subtypes <- my_subtypes[my_subtypes %in% names(clrs$cluster_label[[coi]])]

if (cell_sort == "CD45+") {
cells_to_keep <- colnames(seu_obj)[seu_obj$cluster_label %in% my_subtypes & !(str_detect(seu_obj$cell_id, "CD45N"))]
}

if (cell_sort == "CD45-") {
cells_to_keep <- colnames(seu_obj)[seu_obj$cluster_label %in% my_subtypes & !(str_detect(seu_obj$cell_id, "CD45P"))]
}

# seu_obj_sub <- subset(seu_obj, cells = cells_to_keep)
# seu_obj_sub <- RunUMAP(seu_obj_sub, dims = 1:50, reduction = "harmony", reduction.name = "umapharmony", reduction.key = "umapharmony_")
# seu_obj_sub$cluster_label <- seu_obj$cluster_label[colnames(seu_obj) %in% colnames(seu_obj_sub)]
# write_rds(seu_obj_sub, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_processed_filtered.rds"))
seu_obj_sub <- read_rds(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_processed_filtered.rds"))

plot_data_sub <- as_tibble(FetchData(seu_obj_sub, c(myfeatures, "cluster_label"))) %>% 
  left_join(meta_tbl, by = "sample") %>% 
  mutate(patient_id = str_remove_all(patient_id, "SPECTRUM-OV-"),
         tumor_supersite = ordered(tumor_supersite, levels = names(clrs$tumor_supersite))) %>% 
  mutate(cell_id = colnames(seu_obj_sub),
         cluster_label = ordered(cluster_label, levels = my_subtypes),
         )
  
ggplot(plot_data_sub, aes(umapharmony_1, umapharmony_2, color = cluster_label)) + 
  common_layers_disc +
  ggtitle("Sub cluster") +
  #facet_wrap(~cluster_label) +
  scale_color_manual(values = clrs$cluster_label[[coi]])

ggplot(plot_data_sub, aes(umapharmony_1, umapharmony_2, color = patient_id_short)) + 
  common_layers_disc +
  ggtitle("Patient") +
  # facet_wrap(~therapy) +
  scale_color_manual(values = clrs$patient_id_short)

ggplot(plot_data_sub, aes(umapharmony_1, umapharmony_2, color = tumor_supersite)) + 
  # geom_point(aes(umapharmony_1, umapharmony_2), 
  #            color = "grey90", size = 0.01, 
  #            data = select(plot_data_sub, -tumor_supersite)) +
  common_layers_disc +
  ggtitle("Site") +
  # facet_wrap(~therapy) +
  scale_color_manual(values = clrs$tumor_supersite)

write_tsv(select(plot_data_sub, cell_id, everything(), -umapharmony_1, -umapharmony_2, -contains("RNA_")), paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_embedding.tsv"))

4 sub clustering

big_clusters <- list(
  CD4.T.CXCL13 = "CD4.T.CXCL13", 
  CD4.T.naive = "CD4.T.naive", 
  CD4.T.reg = "CD4.T.reg", 
  CD4.Th17 = "CD4.Th17", 
  CD8.T.CXCL13 = "CD8.T.CXCL13", 
  CD8.T.cytotoxic = "CD8.T.cytotoxic", 
  CD8.T.ISG = "CD8.T.ISG", 
  Cycling.T.NK = "Cycling.T.NK", 
  NK.CD56 = "NK.CD56", 
  NK.cytotoxic = "NK.cytotoxic"
)

cells_to_keep_list <- lapply(big_clusters, function(x) 
colnames(seu_obj_sub)[seu_obj_sub$cluster_label %in% x])

# seu_list <- lapply(cells_to_keep_list, function(x) subset(seu_obj_sub, cells = x))
# 
# preprocess_wrapper <- . %>%
#   FindNeighbors(reduction = "harmony", dims = 1:50) %>%
#   # FindClusters(res = 0.1) %>%
#   # FindClusters(res = 0.2) %>%
#   FindClusters(res = 0.3)
# 
# seu_list <- lapply(seu_list, preprocess_wrapper)
# write_rds(seu_list, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_processed_filtered_sub_clusters.rds"))

seu_list <- read_rds(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_processed_filtered_sub_clusters.rds"))

get_cluster_vector <- function(cluster_label, cluster_res = RNA_snn_res.0.3) {
  cluster_res <- enquo(cluster_res)
  cbind(cell_id = colnames(seu_list[[cluster_label]]), FetchData(seu_list[[cluster_label]], c(as_label(cluster_res)))) %>%
    as_tibble %>%
    mutate(cluster_extended = paste0(cluster_label, "_", !!cluster_res)) %>%
    select(cell_id, cluster_extended) %>% 
    deframe
}

cluster_vector_list <- lapply(big_clusters, get_cluster_vector)

cluster_vector <- unlist(cluster_vector_list, use.names = F) %>% setNames(lapply(cluster_vector_list, names) %>% unlist(use.names = F))

seu_obj$cluster_extended <- cluster_vector[seu_obj$cell_id] %>%
  setNames(seu_obj$cell_id)

seu_obj_sub$cluster_extended <- cluster_vector[seu_obj_sub$cell_id] %>%
  setNames(seu_obj_sub$cell_id)

cluster_extended_uniq <- as.character(na.omit(unique(seu_obj_sub$cluster_extended)))

Idents(seu_obj_sub) <- seu_obj_sub$cluster_extended
Idents(seu_obj) <- seu_obj$cluster_extended

patient_specific_subcluster <- FetchData(seu_obj_sub, c("patient_id", "cluster_extended")) %>% 
  as_tibble() %>% 
  group_by(patient_id, cluster_extended) %>% 
  tally %>% 
  group_by(cluster_extended) %>% 
  mutate(nrel = n/sum(n),
         ntotal = sum(n)) %>% 
  arrange(desc(nrel)) %>% 
  distinct(cluster_extended, .keep_all = T) %>% 
  ungroup

# marker_tbl_extended <- lapply(
#   cluster_extended_uniq,
#   function(x) {
#     marker_tbl <- FindMarkers(seu_obj_sub, ident.1 = x) %>%
#       as_tibble(rownames = "gene") %>%
#       separate(gene, into = c("gene", "drop"), sep = "\\.\\.\\.") %>%
#       select(-drop) %>% 
#       mutate(cluster_extended = x)
#     write_tsv(marker_tbl, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_marker_table_subcluster_", x, ".tsv"))
#     return(marker_tbl)
# }) %>% bind_rows
#  
# write_tsv(marker_tbl_extended, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_marker_table_annotated_extended.tsv"))

marker_tbl_extended <- read_tsv(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_marker_table_annotated_extended.tsv"))
 
cluster_label_extended <- c(
  CD4.T.naive_0 = "CD4.T.naive.centr.mem", 
  CD4.T.naive_1 = "CD4.T.effector.memory", 
  CD4.T.CXCL13_0 = "CD4.T.dysfunc.early", 
  CD4.T.CXCL13_1 = "CD4.T.dysfunc.late.1", 
  CD4.T.CXCL13_2 = "CD4.T.dysfunc.late.2", 
  CD4.T.reg_0 = "CD4.T.reg.2", 
  CD4.T.reg_1 = "CD4.T.reg.1", 
  CD4.T.reg_2 = "CD4.T.reg.ISG", 
  CD4.T.reg_3 = "CD4.T.reg.3", 
  CD4.Th17_0 = "CD4.Th17.1", 
  CD4.Th17_1 = "CD4.Th17.2", 
  CD8.T.cytotoxic_0 = "CD8.T.effector.memory", 
  CD8.T.cytotoxic_1 = "CD8.T.cytotoxic", 
  CD8.T.cytotoxic_2 = "dissociated_3", 
  CD8.T.cytotoxic_3 = "CD8.T.naive.centr.mem", 
  CD8.T.cytotoxic_4 = "gd.T.cell", 
  CD8.T.cytotoxic_5 = "doublet.Plasma.cell_1", 
  CD8.T.CXCL13_0 = "CD8.T.dysfunc.early", 
  CD8.T.CXCL13_1 = "CD8.T.dysfunc.late", 
  CD8.T.CXCL13_2 = "NK.reg.CRTAM", 
  CD8.T.ISG_0 = "CD8.T.ISG", 
  CD8.T.ISG_1 = "CD4.T.ISG", 
  CD8.T.ISG_2 = "dissociated_4", 
  NK.CD56_0 = "NK.reg.IGFBP2",
  NK.CD56_1 = "NK.reg.CD56", 
  NK.CD56_2 = "NK.reg.ISG", 
  NK.CD56_3 = "NK.reg.KRT81.KRT86", 
  NK.CD56_4 = "NK.reg.CCL3",
  NK.cytotoxic_0 = "NK.cytotoxic.SPON2.2", 
  NK.cytotoxic_1 = "NK.cytotoxic.GZMH", 
  NK.cytotoxic_2 = "NK.cytotoxic.SPON2.1",
  Cycling.T.NK_0 = "Cycling.CD8.T.1",
  Cycling.T.NK_1 = "Cycling.CD8.T.2", 
  Cycling.T.NK_2 = "Cycling.CD8.T.3", 
  Cycling.T.NK_3 = "Cycling.CD8.T.4", 
  Cycling.T.NK_4 = "Cycling.NK.3", 
  Cycling.T.NK_5 = "Cycling.CD4.T", 
  Cycling.T.NK_6 = "Cycling.NK.2", 
  Cycling.T.NK_7 = "Cycling.NK.1"
)

names(seu_obj$cluster_label) <- colnames(seu_obj)
seu_obj$cluster_label[is.na(seu_obj$cluster_label)] <- "NA"
Idents(seu_obj) <- seu_obj$cluster_label
seu_obj$cluster_label_sub <- cluster_label_extended[seu_obj$cluster_extended]
names(seu_obj$cluster_label_sub) <- colnames(seu_obj)
seu_obj$cluster_label_sub[is.na(seu_obj$cluster_label_sub)] <- "NA"
Idents(seu_obj) <- seu_obj$cluster_label_sub

names(seu_obj_sub$cluster_label) <- colnames(seu_obj_sub)
seu_obj_sub$cluster_label[is.na(seu_obj_sub$cluster_label)] <- "NA"
Idents(seu_obj_sub) <- seu_obj_sub$cluster_label
seu_obj_sub$cluster_label_sub <- cluster_label_extended[seu_obj_sub$cluster_extended]
names(seu_obj_sub$cluster_label_sub) <- colnames(seu_obj_sub)
seu_obj_sub$cluster_label_sub[is.na(seu_obj_sub$cluster_label_sub)] <- "NA"
Idents(seu_obj_sub) <- seu_obj_sub$cluster_label_sub

marker_sheet_extended <- marker_tbl_extended %>%
  mutate(cluster_label_sub = cluster_label_extended[cluster_extended]) %>% 
  group_by(cluster_label_sub) %>%
  mutate(rank = row_number()) %>% 
  slice(1:50) %>% 
  select(rank, gene, cluster_label_sub) %>% 
  arrange(cluster_label_sub) %>% 
  spread(cluster_label_sub, gene)

cluster_n_tbl_full <- seu_obj$cluster_label_sub %>% 
  table() %>% 
  enframe("cluster_label_sub", "cluster_n") %>% 
  mutate(cluster_nrel = cluster_n/sum(cluster_n)) %>% 
  filter(cluster_label_sub != "NA")

marker_tbl_extended_annotated <- marker_tbl_extended %>% 
  mutate(cluster_label_sub = cluster_label_extended[cluster_extended]) %>%
  select(-cluster_extended) %>% 
  left_join(cluster_n_tbl_full, by = "cluster_label_sub") %>% 
  arrange(cluster_label_sub, desc(avg_logFC))

if (cell_sort == "CD45+") {
cells_to_keep <- colnames(seu_obj_sub)[!(str_detect(seu_obj_sub$cluster_label_sub, "dissociated|doublet")) & !(str_detect(seu_obj_sub$cell_id, "CD45N"))]
}

if (cell_sort == "CD45-") {
cells_to_keep <- colnames(seu_obj_sub)[!(str_detect(seu_obj_sub$cluster_label_sub, "dissociated|doublet")) & !(str_detect(seu_obj_sub$cell_id, "CD45P"))]
}


# seu_obj_sub_sub <- subset(seu_obj_sub, cells = cells_to_keep)
# seu_obj_sub_sub <- RunUMAP(seu_obj_sub_sub, dims = 1:50, reduction = "harmony", reduction.name = "umapharmony", reduction.key = "umapharmony_")
# write_rds(seu_obj_sub_sub, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_processed_filtered_sub.rds"))

seu_obj_sub_sub <- read_rds(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_processed_filtered_sub.rds"))

seu_obj_sub_sub$cluster_label_sub <- setNames(cluster_label_extended[seu_obj_sub_sub$cluster_extended], seu_obj_sub_sub$cell_id)

write_rds(seu_obj_sub_sub, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_processed_filtered_sub.rds"))

formattable::formattable(marker_sheet_extended)
rank CD4.T.dysfunc.early CD4.T.dysfunc.late.1 CD4.T.dysfunc.late.2 CD4.T.effector.memory CD4.T.ISG CD4.T.naive.centr.mem CD4.T.reg.1 CD4.T.reg.2 CD4.T.reg.3 CD4.T.reg.ISG CD4.Th17.1 CD4.Th17.2 CD8.T.cytotoxic CD8.T.dysfunc.early CD8.T.dysfunc.late CD8.T.effector.memory CD8.T.ISG CD8.T.naive.centr.mem Cycling.CD4.T Cycling.CD8.T.1 Cycling.CD8.T.2 Cycling.CD8.T.3 Cycling.CD8.T.4 Cycling.NK.1 Cycling.NK.2 Cycling.NK.3 dissociated_3 dissociated_4 doublet.Plasma.cell_1 gd.T.cell NK.cytotoxic.GZMH NK.cytotoxic.SPON2.1 NK.cytotoxic.SPON2.2 NK.reg.CCL3 NK.reg.CD56 NK.reg.CRTAM NK.reg.IGFBP2 NK.reg.ISG NK.reg.KRT81.KRT86
1 CXCL13 CXCL13 CXCL13 IL7R IFIT3 CCR7 RTKN2 TNFRSF4 TNFRSF4 ISG15 KLRB1 LST1 GZMH CXCL13 CXCL13 GZMK IFIT3 CD8B STMN1 STMN1 STMN1 CXCL13 CXCL13 SPON2 CENPF CENPF GZMK XIST IGHG3 TRGC2 FGFBP2 SPON2 FGFBP2 CCL3 GNLY XCL2 TYROBP ISG15 GNLY
2 MAF NMB PDCD1 CD40LG ISG15 RPS8 FOXP3 IL2RA TNFRSF18 SAT1 IL7R IL4I1 CD8A CCL4L2 KRT86 CD8A ISG15 IL7R MKI67 TUBA1B MKI67 STMN1 HLA-DRA FCGR3A ASPM MKI67 MALAT1 MALAT1 IGHG1 XCL1 GZMH FGFBP2 FCGR3A AREG AREG GZMB KLRC1 GNLY KLRC1
3 IL6ST CPM CTLA4 ANXA1 MX2 RPS6 LTB FOXP3 TNFRSF9 HERC5 IL4I1 TNFSF13B CD8B HAVCR2 TNFRSF9 CD8B IFIT1 CD8A TUBA1B HIST1H4C TOP2A TYMS HLA-DRB1 IGFBP7 MKI67 STMN1 PTPRC NEAT1 IGLC3 ZNF683 CX3CR1 FCGR3A SPON2 CCL4L2 KLRC1 CRTAM FCER1G IFIT3 KRT81
4 NR3C1 NR3C1 MAF GPR183 IFIT2 RPS3A UGP2 BATF IL2RA MX1 LTB PCDH9 CCL5 MIR155HG GZMB KLRG1 MX1 YBX3 HMGB2 MKI67 ASPM TUBA1B GZMB KLRF1 TOP2A TOP2A EDF1 MT-ND3 IGLC2 KLRC2 NKG7 CCL3 KLRF1 XCL2 ITGAX XCL1 TRDC FCER1G FCER1G
5 FKBP5 IGFL2 ZBED2 KLF2 IFI44L SELL TBC1D4 TNFRSF18 BATF TNFRSF4 CCR6 SPINK2 HLA-DRB1 PTMS RBPJ ITM2C IFIT2 ANXA1 TOP2A TUBB TYMS PCNA IFNG CX3CR1 KNL1 ASPM RPS29 MT-ND1 IGHG4 LINC02446 TRGC2 PRF1 PRF1 XCL1 FCER1G FABP5 GNLY MX1 AREG
6 CD40LG FKBP5 RBPJ S100A4 IFIT1 EEF1B2 IL2RA LTB CTLA4 SPATS2L TNFSF13B KIT GZMA CD8A SRGAP3 CCL4 RSAD2 RPS12 TUBB TYMS HELLS PCLAF FABP5 AKR1C3 NUSAP1 TUBA1B RPL36AL MT-CO3 IGKC HOPX KLF2 KLRF1 GNLY FOS XCL1 TNFRSF9 AREG IFIT1 TNFRSF18
7 CD4 ITM2A DUSP4 FOS RSAD2 RPS12 PMAIP1 CTLA4 GADD45A IFI6 LST1 SCN1B GZMK LAG3 FAM3C DUSP2 IFI6 TOB1 CENPF TOP2A CLSPN DUT COTL1 GINS2 CENPE HMGB2 RPL13A MT-ND2 JCHAIN XCL2 KLRG1 PTGDS KLRD1 IER2 TYROBP HSP90AB1 IGFBP2 TYROBP KRT86
8 AC004585.1 ICA1 CD4 LTB MX1 RPL34 IKZF2 GADD45A TNFRSF1B IL2RA AQP3 LTC4S HLA-DPB1 GZMB TIGIT PPP1R14B ISG20 FOS HMGN2 HMGB2 PCNA FABP5 GAPDH FGFBP2 TPX2 GNLY ATP5MC2 MT-CYB IGHGP CTSW PLEK CX3CR1 NKG7 CCL4 TRDC PKM KRT81 IFI27 TRDC
9 CORO1B IL6ST CD40LG TPT1 IL7R RPL32 IL32 RTKN2 ICOS CTLA4 CCL20 IL1R1 PTMS IFNG SPRY1 TRGC2 MX2 VIM PTTG1 NUSAP1 UBE2C CLSPN NASP MLC1 DLGAP5 TUBB CALM1 MT-ATP6 IGHG2 CD63 PRF1 GZMB CLIC3 EGR1 XCL2 GAPDH KLRD1 IFITM3 TYROBP
10 ITM2A TSHZ2 MIR155HG NFKBIA IFI6 KLF2 CTLA4 SOX4 PIM3 OAS1 CTSH LINC00299 ZNF683 CCL3 IFNG CST7 SAMD9L JUNB IL2RA CENPF PCLAF HELLS LAG3 S1PR5 GTSE1 NUSAP1 SRP14 MT-CO1 CD8A CD7 ZEB2 IGFBP7 TYROBP NFKBIA KLRD1 NME1 CLIC3 KLRC1 XCL1
11 TSHZ2 TNFRSF4 PRDM1 TIMP1 TNFSF10 RPL9 TIGIT TNFRSF1B FOXP3 IFIT3 RORA GSN HLA-DPA1 CXCR6 MIR155HG DTHD1 IFI44L GZMK PCLAF PCLAF GTSE1 NASP MCM5 TYMS KIF14 UBE2C COMMD6 PTPRC MZB1 TRGC1 PRSS23 PLAC8 GZMB CCL3L1 CMC1 MIR155HG KLRB1 AREG CTSW
12 PDCD1 AHI1 TNFRSF18 S100A10 STAT1 IL7R ARID5B TIGIT IL1R2 IFIT1 TNFRSF25 PLAT CXCR6 HLA-DRB1 CD8A CCL4L2 HERC5 CCR7 NUSAP1 HIST1H1B ORC6 MCM3 PTMS E2F1 HJURP IGFBP2 RPS25 MT-CO2 CCL5 NCR3 KLRD1 GNLY PLAC8 DUSP2 CTSW SRM CEBPD IFI44L ZNF683
13 TOX2 AC004585.1 ITM2A JUN HERC5 RPS13 SAT1 TBC1D4 LTB FOXP3 SLC4A10 CD300LF HLA-DRA CCL5 LAG3 CRTAM OAS1 ZFP36L2 ASPM HMGN2 CDT1 MKI67 MIR155HG MCM7 POLQ TPX2 PFDN5 H3F3B IGHA1 CD52 FCGR3A FCER1G CX3CR1 BTG2 NCAM1 NPW KRT86 RSAD2 XCL2
14 CH25H LIMS1 RGS1 S100A11 ISG20 RPL21 BATF LINC01943 TBC1D4 TYMP DPP4 AFF3 ITM2C RBPJ PHLDA1 GZMM MT2A S100A6 H2AFZ ASPM BRCA1 MCM5 HAVCR2 MYOM2 BUB1B KLRC1 RPS15A CCNI CCL4 TRDC MYBL1 NKG7 FCER1G NR4A1 HOPX KDM6B TXK CD38 HOPX
15 PTPN13 BIRC3 LINC01943 VIM OAS1 TCF7 CORO1B SAT1 PKM LTB MYBL1 DLL1 COTL1 GZMH ENTPD1 CMC1 STAT1 MT-ND1 TYMS H2AFZ ZWINT MCM7 CD8A MCM4 HMMR HMGN2 HINT1 GMFG GZMK KLRD1 ADGRG1 AKR1C3 PTGDS FOSB LINC00996 TSPAN17 XCL1 IFIT2 RGS16
16 CTSW PGM2L1 TNFRSF4 CTSL OAS3 MAL MIR4435-2HG TNFRSF9 BIRC3 RTKN2 S100A4 LIF HLA-DQA1 HLA-DRA HAVCR2 LYAR TNFSF10 KLF2 CEP55 RRM2 ATAD2 MCM4 HELLS CHEK1 KIF11 SMC4 RPS27 COX5B CD8B KIR3DL2 EFHD2 AREG PLEK CD69 KLRF1 ZBED2 IL2RB IFI6 CD7
17 CCL5 CD200 ADAM19 RPL34 EIF2AK2 RPL39 TNFRSF4 UGP2 REL EPSTI1 CD40LG IL23R IFNG PHLDA1 CCL3 CCL5 OASL LYAR DUT SMC4 RAD51AP1 TUBB MCM6 CLSPN CIT FCER1G FTH1 CALM2 RGS1 TRG-AS1 FLNA CCL4 EFHD2 TYROBP IFITM3 GOLIM4 MATK TRDC LGALS1
18 NKG7 SESN3 CCDC50 EEF1A1 IFI44 RPL22 AC005224.3 IKZF2 ENTPD1 STAT1 TPT1 RORC CCL4 PDCD1 CTLA4 EOMES PLSCR1 CXCR4 PCNA PCNA NUSAP1 GAPDH ACTB STMN1 SGO2 TYMS HSPA8 PSMB1 GZMA CCL5 GNLY KLRD1 CST7 NR4A2 TXK PAICS NKG7 TNFSF10 CD63
19 KLRD1 RNF19A TSHZ2 RPLP0 GBP1 RPL7 PBXIP1 PMAIP1 TYMP ICOS RPLP1 TGM2 APOBEC3G HLA-DPA1 JAML SH2D1A GBP1 NELL2 TNFRSF4 DUT DLGAP5 GINS2 MCM4 CHST2 C21orf58 GSTP1 DBI TPM3 HLA-DPA1 TMSB4X CST7 PLEK CCL3 GADD45B IL2RB MRTO4 XCL2 ISG20 KLRC2
20 ST8SIA1 CHN1 CD2 RPS18 CCR7 RPL5 STAM IL32 SAT1 ISG20 EEF1A1 TNFSF11 SH3BGRL3 CD8B CCND2 YBX3 SAMD9 S100A10 TPX2 ATAD2 CDK1 HMGB2 STMN1 CTBP2 KIF23 KRT81 FAU SOD1 CXCR6 LINC01871 HOPX CLIC3 IGFBP7 KLRD1 KLRC2 GEM CD63 GZMB TMIGD2
21 LIMS1 GNG4 ETV7 RPS12 EPSTI1 RPS18 GADD45A ICOS IL1R1 MX2 B3GALT2 ZG16B CD52 HLA-DPB1 PTMS GZMA EIF2AK2 MT-ND2 FOXP3 CLSPN BIRC5 RANBP1 PCNA FCER1G AC007240.1 TRDC PARK7 ATP5MPL IKZF3 GNLY C12orf75 EFHD2 ADGRG1 TRDC MATK EGR2 TMIGD2 SAMD9L SRGAP3
22 SMCO4 MAF MYO7A RPL39 MT2A RPL37 BTG3 LAYN SDC4 TIGIT RORC ENPP1 CD3D VCAM1 LINC01871 TC2N LAG3 TNFRSF18 ENO1 UBE2C TUBA1B HLA-DRA CLSPN PLAC8 KIF15 H2AFZ RPLP2 AES LTB ZFP36L2 GZMB MYOM2 HOPX KLRC1 SELL CD72 CMC1 BST2 MATK
23 CD8A BATF TOX2 JUNB SAMD9L RPL13 BIRC3 LINC02099 LINC01943 RSAD2 LTK PTGER3 LAG3 LINC01871 GOLIM4 CXCR4 IFI35 LMNA BATF HMGB1 NASP MCM6 GINS2 RAMP1 NCAPG TUBB4B C9orf16 GNG5 APOC1 HCST C1orf21 TYROBP CD247 FCER1G CD7 METTL1 CTSW CLIC3 CXXC5
24 CTLA4 TOX2 CD200 TOB1 XAF1 RPS5 HPGD CTSC NAMPT TNFRSF1B CFH IRAK3 C12orf75 JAML GAPDH PPP2R5C EPSTI1 TSC22D3 CTLA4 HIST1H1D MCM4 DEK CDCA7 ORC6 BRCA2 AREG RPL23A SAP18 RPL22 IKZF2 LITAF ADGRG1 AREG NFKBID GZMB GFOD1 CD7 PLSCR1 GSTP1
25 CCL4 CH25H BHLHE40-AS1 AQP3 TMEM123 RPL3 SELL DUSP4 LTA UGP2 NCR3 CA2 CD74 TNIP3 TNS3 A2M-AS1 OAS3 FOSB GAPDH TPX2 GINS2 H2AFZ CHEK1 MCM2 ECT2 CENPE RPL28 CIB1 PABPC1 FCRL6 PLAC8 CST7 ZEB2 ZFP36 SRGAP3 RGS16 CCL3 OAS1 IL2RB
26 FYB1 THADA PTPN13 RPL11 DDX58 EEF1A1 ICOS CD27 DNPH1 IFI44L FAM241A KIAA1211L CD3G TNFSF4 LAYN GZMH XAF1 ZFP36 HIST1H1B H2AFX CENPU HMGN2 MCM2 TTC38 CKAP2L UBE2S COX7C COPE HLA-DPB1 IL32 S100A4 HOPX AKR1C3 CEBPD CEBPD PARK7 LAT2 IFI35 LAT2
27 RUNX3 CD84 PDE7B TNFAIP3 CMPK2 RPL11 AC133644.2 ARID5B TIGIT SAMD9 CEBPD SLC4A10 THEMIS CD63 NBL1 PDCD4 USP18 PPP1R14B SMC4 SMC2 FEN1 SMC2 TYMS GZMB ARHGAP11B CKS1B CHCHD2 TRMT112 TPT1 NUCB2 SPON2 CD247 KLRB1 MAFF CD63 NDFIP2 GSTP1 KLRD1 CSF1
28 RNF19A BTLA NKG7 TNFRSF25 CD40LG RPS23 NABP1 CORO1B TFRC TNFRSF18 RPLP0 COL4A4 LINC02446 HLA-DQA1 MYO1E HLA-DPB1 LY6E DUSP1 SMC2 TMPO TUBB DNMT1 CDT1 PCNA ANLN TNFRSF18 RPL4 PSMB9 CD74 CD9 FCRL6 PRSS23 PTPN12 IGFBP2 NKG7 IFNG HOPX EPSTI1 SPRY1
29 PRF1 CORO1B IFNG FXYD5 SAMD9 RPL23 CD27 ENTPD1 F5 SAMD9L RPS12 NRP1 HLA-DQB1 ID2 RGS2 F2R CMPK2 CD55 BIRC5 CKS1B CHEK1 TPI1 MALAT1 PLEK ESPL1 HIST1H1B RPS24 BRK1 GZMH CD8A ADRB2 ABHD17A PRSS23 KDM6B TNFRSF18 COL6A3 CCL5 IRF7 CAPG
30 RGS1 METTL8 NR3C1 RPL3 LY6E TPT1 IKZF4 MAGEH1 PMAIP1 TBC1D4 ERN1 B3GALT5 LGALS1 FABP5 MYO7A SLAMF7 NT5C3A CD69 PKM UBE2S KNL1 ATAD2 HLA-DQA1 MCM5 NUF2 KPNA2 COX8A RBX1 IGHM KLRC4 ARL4C CTSW C1orf21 IER3 IFITM2 CXCL13 KLRF1 CEBPD ITGA1
31 GZMB CCDC50 LIMS1 RPL9 GPR183 RPLP0 PHTF2 S100A4 CCND2 BATF KIF5C PDZK1 HLA-DRB5 LINC02446 AKAP5 TNFSF9 IFI44 GPR183 RRM2 KNL1 SMC2 BRCA1 MYO7A HES6 KIF4A DLGAP5 RPS16 HSPA8 RPL8 MATK S1PR5 ZEB2 TXK CD160 CD300A TNFRSF18 SH2D1B XCL2 CLIC3
32 IL6R SMCO4 RNF19A RPL32 USP18 RPL30 PHACTR2 BIRC3 DUSP4 USP18 RPL13 EMID1 S100A6 TNFRSF9 DGKH STK17A HELZ2 TYROBP PMAIP1 TUBB4B MCM5 HNRNPAB MKI67 PRF1 KIFC1 CKS2 RPL35A S100A11 CD3E CRTAM KLF3 CHST2 S1PR5 REL GSTP1 RGCC PRF1 HERC5 LDLRAD4
33 CST7 CD4 RDH10 RPS3A DDX60L RPS28 RGS1 CD4 LINC02099 COX5A SPOCK2 HOXA5 FABP5 TIGIT SNX9 PKM IRF7 FCER1G TK1 HIST1H1C KIF11 SLBP HLA-DPA1 FEN1 BUB1 XCL2 RPL7 VAMP8 HLA-DQB1 S100A6 PATL2 PTPN12 CEBPD GNLY KRT81 CD82 ITGA1 NT5C3A KIR2DL4
34 SLC9A9 ARMH1 GAPDH RPS6 TRIM22 LEF1 DUSP16 MAF CCL20 SOX4 JAML LINGO4 ANXA5 ITGAE DUSP4 COX5A CD38 CCL3 JPT1 EZH2 TPX2 SNRPB SH3BGRL3 CDT1 KIF20B SMC2 RPS7 POMP WNK1 CAPG ANXA1 MYBL1 TRDC KLRF1 SPTSSB MALAT1 SRGAP3 IL2RB IGFBP2
35 BTLA PHACTR2 SNX9 CD4 HELZ2 RPS20 TRAC SPOCK2 SH2D2A IKZF2 TIGIT NECTIN2 CCL4L2 LYST LINC01480 EIF2AK2 TRIM22 CTSC ANP32B DEK HIST1H1B GZMB MCM3 ADGRG1 FBXO43 TMPO RPS11 EDF1 CD4 IFITM2 TGFBR3 TTC38 CTSW CLIC3 IGFBP2 RBPJ MCTP2 LY6E TXK
36 TP53INP1 TNFRSF25 IL6ST RPL36 ANXA1 RPL38 CARD16 IL1R1 ZBTB32 PMAIP1 RPL34 HOXA7 TMSB4X APOBEC3C ITGAE ACP5 DDX58 RGCC TMPO CDK1 GMNN COTL1 PTTG1 CCDC28B NCAPG2 CDK1 ATP5MG NDUFS5 HLA-DRA MAF CTSW S1PR5 GZMH TXK CAPN12 RANBP1 MAPK1 XAF1 LINC01871
37 TNFSF8 ARID5B CORO1B RPL30 NT5C3A RPS16 LINC01943 GLRX ID3 CTSC S100A6 PRAM1 LGALS3 CD27 RGS1 ADGRG1 PARP14 AC020916.1 CENPE BIRC5 MCM7 DNAJC9 MCM7 GNLY PRC1 CTSW MYL12B UQCR11 RPS13 CD4 S100A6 ITGB2 TTC38 MATK CXXC5 EIF5A CXXC5 XCL1 CD9
38 GZMH CTLA4 ITGA2 RPS8 PLSCR1 RPS4X S100A4 CARD16 SOX4 NT5C3A RPL32 B4GALNT1 ITGA1 HLA-DMA CD63 IFI35 CD8A CTLA4 UBE2C H2AFV RRM1 HLA-DRB1 RBPJ ASCL2 DIAPH3 KNL1 RPL39 LDHB CD3G GLUL TTC16 XBP1 ABHD17A KLRB1 IL18 DCTPP1 KLRC2 PRF1 SLC16A3
39 GNLY ZBED2 SMCO4 RPS27A IFIH1 RPL36 TNFRSF1B TYMP PTP4A3 TNFRSF9 HLA-DPB1 LPAR1 TRGC2 CRTAM SNAP47 CSF1 DDX60 CD4 CLSPN TK1 WDR76 CENPX CD8B PTPN12 SMC4 KRT86 UBB HNRNPA1 RPL29 PRMT9 A2M-AS1 CD160 PTGDR NKG7 ID2 ENO1 SAMD3 SAMD9 CCL5
40 METTL8 TNFRSF18 SAMSN1 RPSA IFI35 RPLP2 MAGEH1 MIR4435-2HG MIR155HG IFI44 TLE1 SLCO2A1 APOBEC3C NKG7 CLNK NDFIP2 PPM1K DNAJB1 CKS2 CENPM RRM2 IFNG PDCD1 PCLAF RPL10 PCLAF RPS21 ANP32B TNFRSF18 SCML4 TTC38 LAIR2 MYBL1 EGR2 PTGDR TPI1 LINC00996 MX2 NCR3
41 ADAM19 RILPL2 SRGN RPL13 IRF7 RPS27 CD4 DNPH1 ARID5B ARID5B TMIGD2 CD33 SIT1 CTLA4 TNIP3 LGALS3 IFIH1 MYBL1 ACTB ZWINT CENPF CHEK1 BRCA1 CDC45 RPS27 TYROBP RPL13 ABRACL LARP7 LYAR ASCL2 KLRB1 CMC1 IL2RB KLRB1 CAMK1 CD247 OASL PRMT9
42 CD2 AGFG1 ICOS EEF1B2 DDX60 RPS2 SPOCK2 PHACTR2 CD27 LINC01943 FKBP11 APOL4 JAML RGS1 ADAMTS6 TMEM123 RNF213 CD38 CALM3 HIST1H1E UBE2T CCL3 ITGA2 GMNN BRCA1 XCL1 ANP32B TRAPPC1 DUSP23 CTSA VCL TXK CHST2 MAP3K8 SH2D1B TNIP3 PTGDR USP18 CLNK
43 SPOCK2 CYSLTR1 KLRD1 RPL38 OAS2 RPL12 CCNG2 IL1R2 ENO1 GBP5 HLA-DRB1 RASSF8 LINC01871 ACP5 BCL2L11 IL6ST DDX60L TNFRSF4 ATAD2 FAM111B CDCA7 CDT1 CD63 CENPU RPL41 HMMR COX6C ATP5F1D DUSP2 IFITM3 LILRB1 PTGDR CCL4 IRF8 C1orf162 PIM3 NCAM1 EIF2AK2 CTSA
44 GZMA CD40LG CD109 RPLP1 PNPT1 SNHG8 FCMR PHLDA1 IKZF2 OAS3 RPS27A TLE1 CD2 SNAP47 NDFIP2 CD247 OAS2 TNFAIP3 ANP32E KIF11 CDKN3 RAN CCL5 C1orf21 RPL30 TMIGD2 RPS28 SNRPB SELL PTGDR FGR CD300A SH2D1B KRT81 PLAC8 NAMPT STARD3NL ERICH3 KLRB1
45 CD8B BHLHE40-AS1 CD3D RPL5 PPM1K RPS27A TLK1 STAM CD4 CD4 KLRD1 AL034397.3 CST7 ITM2A RAB27A TPI1 PARP9 HLA-DRA ICOS FABP5 SGO2 PKM CCL3 FGR RPL13 CDKN3 RPL36 SEPT7 KLF2 SPRY1 PLEKHG3 FGR CEP78 SRGN KIR2DL4 TNFRSF1B CAPN12 CCL3 CTSD
46 MATK IGFBP4 SLAMF1 RPS16 PARP9 RPL19 SMAP2 SLAMF1 LAIR2 RNF213 RUNX2 STAC CCR5 SAMSN1 TNFRSF1B ANKRD28 PNPT1 JUN CKS1B DLGAP5 USP1 ENO1 HLA-DMA LYN RPS27A H2AFX LDHB ZFAS1 CXCR4 CTLA4 COL6A2 GSAP XBP1 LAT2 MCTP2 ZNF593 TNFRSF18 IFI44 NDFIP2
47 CD84 ELMO1 RPL23A CD52 LGALS3BP PABPC1 LAYN LAIR2 SYNGR2 SPOCK2 RPL10 BCAS1 CRTAM CCND2 ID2 CD226 GZMK NCF1 PHF19 GTSE1 FANCI SNRPD1 GZMH KIR2DL3 RPS12 CD63 SUMO2 UQCRH DNAJC7 DDAH2 CES1 CTBP2 SAMD3 ZBTB16 ZNF683 HSPA5 AOAH CD7 ADGRG3
48 ARID5B MS4A6A AHI1 LMNA LTB RPS14 PIM2 F5 SERPINB9 ENTPD1 HLA-DPA1 CSF2 GRAP2 APOBEC3G LYST TMEM173 BST2 TRDC TNFRSF18 TMEM106C MT-ND4 NUDT1 APOBEC3G ITGAM PRR11 CCNB1 RPL14 SEC61G TRBC2 GUK1 LINC02384 GK5 MATK GRASP MMP25-AS1 FAM3C GZMA OAS3 CD247
49 CD63 ST8SIA1 ARID5B RPL37 LAMP3 RPS21 SAMHD1 TRAC NR4A1 CD27 RPS13 TOX2 IDH2 CD74 AHI1 TPM4 SP110 HAVCR2 UBE2S ANP32B ATAD5 PGAM1 ZNF683 MCM10 RPS14 GTSE1 RPS19 PPP1CA SERPINB6 CD28 COTL1 C1orf21 ITGB2 NR4A3 MAFF SLC7A5 APOBEC3G KLRC2 GZMB
50 HOPX MAGEH1 TYMP RPS29 GZMA RPL35A ICA1 PHTF2 VDR LY6E GZMB ALDOC GZMM DUSP4 ZEB2 USP18 C19orf66 MYADM RANBP1 HNRNPAB CENPE ATAD5 ANXA5 CD300A TPT1 BIRC5 ATP5F1E RHOA TYROBP CCR7 CD27 MTSS1 CD300A EGR3 XBP1 PGAM1 METRNL WARS CKLF
write_tsv(marker_sheet_extended, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_marker_sheet_full.tsv"))

write_tsv(marker_sheet_extended, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/supplementary_tables/", coi, "_marker_sheet_full.tsv"))

write_tsv(marker_tbl_extended_annotated, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_marker_table_annotated_full.tsv"))

write_tsv(marker_tbl_extended_annotated, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/supplementary_tables/", coi, "_marker_table_annotated_full.tsv"))
plot_data_sub_sub <- as_tibble(FetchData(seu_obj_sub_sub, c(myfeatures, "cluster_label_sub"))) %>% 
  left_join(meta_tbl, by = "sample") %>% 
  mutate(cell_id = colnames(seu_obj_sub_sub),
         cluster_label_sub = ordered(cluster_label_sub, levels = names(clrs$cluster_label_sub[[coi]])),
         ) %>%
  mutate(patient_id = str_remove_all(patient_id, "SPECTRUM-OV-"),
         tumor_supersite = ordered(tumor_supersite, levels = names(clrs$tumor_supersite)))
  
if (cell_sort == "CD45+") {
plot_data_sub_sub <- filter(plot_data_sub_sub, sort_parameters != "singlet, live, CD45-", !is.na(tumor_supersite))
}

if (cell_sort == "CD45-") {
plot_data_sub_sub <- filter(plot_data_sub_sub, sort_parameters != "singlet, live, CD45+", !is.na(tumor_supersite))
}

ggplot(plot_data_sub_sub, aes(umapharmony_1, umapharmony_2, color = cluster_label_sub)) + 
  common_layers_disc +
  # facet_wrap(~cluster_label_sub) +
  scale_color_manual(values = clrs$cluster_label_sub[[coi]]) +
  ggtitle("Sub cluster") 

4.1 add module scores and pathway scores

signature_modules_cd8 <- read_excel("_data/small/signatures/SPECTRUM v5 sub cluster markers.xlsx", sheet = 2, skip = 1, range = "M2:P100") %>% 
  gather(module, gene) %>% 
  na.omit() %>% 
  group_by(module) %>% 
  do(gene = c(.$gene)) %>% 
  {setNames(.$gene, .$module)} %>% 
  setNames(paste0(names(.), ".module"))

signature_modules_cin <- yaml::read_yaml("_data/small/signatures/cgas_sting.yaml") %>% 
  .[lapply(., length) > 1] %>% 
  setNames(paste0(names(.), ".module"))

## compute expression module scores
for (i in 1:length(signature_modules_cd8)) {
  seu_obj_sub_sub <- AddModuleScore(seu_obj_sub_sub, features = signature_modules_cd8[i], name = names(signature_modules_cd8)[i])
  seu_obj_sub_sub[[names(signature_modules_cd8)[i]]] <- seu_obj_sub_sub[[paste0(names(signature_modules_cd8)[i], "1")]]
  seu_obj_sub_sub[[paste0(names(signature_modules_cd8)[i], "1")]] <- NULL
  print(paste(names(signature_modules_cd8)[i], "DONE"))
}
## [1] "CD8.Cytotoxic.module DONE"
## [1] "CD8.Dysfunctional.module DONE"
## [1] "CD8.Naive.module DONE"
## [1] "CD8.Predysfunctional.module DONE"
for (i in 1:length(signature_modules_cin)) {
  seu_obj_sub_sub <- AddModuleScore(seu_obj_sub_sub, features = signature_modules_cin[i], name = names(signature_modules_cin)[i])
  seu_obj_sub_sub[[names(signature_modules_cin)[i]]] <- seu_obj_sub_sub[[paste0(names(signature_modules_cin)[i], "1")]]
  seu_obj_sub_sub[[paste0(names(signature_modules_cin)[i], "1")]] <- NULL
  print(paste(names(signature_modules_cin)[i], "DONE"))
}
## [1] "CIN.module DONE"
## [1] "CIN.responsive.Noncanonical.NFKB.targets.up.module DONE"
## [1] "CIN.responsive.Noncanonical.NFKB.targets.down.module DONE"
## [1] "Noncanonical.NFKB.regulators.up.module DONE"
## [1] "Noncanonical.NFKB.regulators.down.module DONE"
## [1] "Canonical.NFKB.regulators.module DONE"
## [1] "Interferon.regulators.module DONE"
## [1] "EMT.module DONE"
## [1] "Inflammation.module DONE"
## [1] "Motility.and.migration.module DONE"
## [1] "ISG.module DONE"
## [1] "SASP.module DONE"
## [1] "Myofibroblast.module DONE"
## [1] "Fibroblast.module DONE"
## [1] "Mesenchymal.stem.cell.module DONE"
## [1] "Pericytes.module DONE"
## [1] "Smooth.muscle.cell.module DONE"
## [1] "ER.stress.module DONE"
## [1] "IFNg.signaling.module DONE"
## compute progeny scores
progeny_list <- seu_obj_sub_sub@assays$RNA@data[VariableFeatures(seu_obj_sub_sub),] %>%
  as.matrix %>%
  progeny %>%
  as.data.frame %>%
  as.list

names(progeny_list) <- make.names(paste0(names(progeny_list), ".pathway"))

for (i in 1:length(progeny_list)) {
  seu_obj_sub_sub <- AddMetaData(seu_obj_sub_sub,
                                 metadata = progeny_list[[i]],
                                 col.name = names(progeny_list)[i])
}

write_rds(seu_obj_sub_sub, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_processed_filtered_annotated.rds"))

seu_obj_sub_sub <- read_rds(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_processed_filtered_annotated.rds"))

seu_obj_sub_sub$cluster_label_sub <- setNames(cluster_label_extended[seu_obj_sub_sub$cluster_extended], seu_obj_sub_sub$cell_id)

write_rds(seu_obj_sub_sub, paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_processed_filtered_annotated.rds"))

4.2 marker heatmap

marker_top_tbl <- marker_sheet_extended[,-1] %>% 
  slice(1:10) %>% 
  as.list %>% 
  .[!str_detect(names(.), "doublet|dissociated")] %>% 
  enframe("cluster_label_x", "gene") %>% 
  unnest(gene)

plot_data_markers <- as_tibble(FetchData(seu_obj_sub_sub, c("cluster_label_sub", myfeatures, unique(marker_top_tbl$gene)))) %>% 
  gather(gene, value, -c(1:(length(myfeatures)+1))) %>% 
  left_join(meta_tbl, by = "sample") %>% 
  mutate(patient_id = str_remove_all(patient_id, "SPECTRUM-OV-"),
         tumor_supersite = ordered(tumor_supersite, levels = names(clrs$tumor_supersite))) %>% 
  mutate(cluster_label_sub = ordered(cluster_label_sub, levels = names(clrs$cluster_label_sub[[coi]]))) %>% 
  group_by(cluster_label_sub, gene) %>% 
  summarise(value = mean(value, na.rm = T)) %>% 
  group_by(gene) %>% 
  mutate(value = scales::rescale(value)) %>% 
  left_join(marker_top_tbl, by = "gene") %>% 
  mutate(cluster_label_x = ordered(cluster_label_x, levels = rev(names(clrs$cluster_label_sub[[coi]]))))

ggplot(plot_data_markers) +
  geom_tile(aes(gene, cluster_label_sub, fill = value)) +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5)) +
  facet_wrap(~cluster_label_x, scales = "free") +
  scale_fill_gradientn(colors = viridis(9)) +
  labs(fill = "Scaled\nexpression") +
  theme(aspect.ratio = 1,
        axis.line = element_blank(),
        axis.ticks = element_blank(),
        axis.title = element_blank())

# ggsave(paste0("_fig/002_marker_heatmap_", coi, ".pdf"), width = nrow(marker_top_tbl)/6, height = 5)

4.3 composition

4.3.1 per site

comp_site_tbl <- plot_data_sub_sub %>%
  filter(!is.na(tumor_supersite)) %>% 
  group_by(cluster_label_sub, tumor_supersite) %>%
  tally %>%
  group_by(tumor_supersite) %>%
  mutate(nrel = n/sum(n)*100) %>%
  ungroup

pnrel_site <- ggplot(comp_site_tbl) +
  geom_bar(aes(tumor_supersite, nrel, fill = cluster_label_sub),
           stat = "identity", position = position_stack()) +
  scale_y_continuous(expand = c(0, 0)) +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5)) +
  labs(fill = "Cluster", y = "Fraction [%]", x = "") +
  scale_fill_manual(values = clrs$cluster_label_sub[[coi]])

pnabs_site <- ggplot(comp_site_tbl) +
  geom_bar(aes(tumor_supersite, n, fill = cluster_label_sub),
           stat = "identity", position = position_stack()) +
  scale_y_continuous(expand = c(0, 0)) +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5)) +
  labs(fill = "Cluster", y = "# cells", x = "") +
  scale_fill_manual(values = clrs$cluster_label_sub[[coi]])

plot_grid(pnabs_site, pnrel_site, ncol = 2, align = "h")

# ggsave(paste0("_fig/02_deep_dive_", coi, "_comp_site.pdf"), width = 8, height = 4)

4.3.2 per sample

comp_tbl_sample_sort <- plot_data_sub_sub %>% 
  group_by(tumor_subsite, tumor_supersite, patient_id, consensus_signature, therapy, cluster_label_sub) %>% 
  tally %>% 
  group_by(tumor_subsite, tumor_supersite, patient_id, consensus_signature, therapy) %>% 
  mutate(nrel = n/sum(n)*100,
         nsum = sum(n),
         log10n = log10(n),
         label_supersite = "Site",
         label_mutsig = "Signature",
         label_therapy = "Rx") %>% 
  ungroup %>% 
  arrange(desc(therapy), tumor_supersite) %>% 
  mutate(tumor_subsite_rx = paste0(tumor_subsite, "_", therapy)) %>% 
  mutate(tumor_subsite = ordered(tumor_subsite, levels = unique(tumor_subsite)),
         tumor_subsite_rx = ordered(tumor_subsite_rx, levels = unique(tumor_subsite_rx))) %>% 
  arrange(patient_id) %>% 
  mutate(label_patient_id = ifelse(as.logical(as.numeric(fct_inorder(as.character(patient_id)))%%2), "Patient1", "Patient2"))

sample_id_x_tbl <- plot_data_sub %>% 
  mutate(sort_short_x = cell_sort) %>% 
  distinct(patient_id, sort_short_x, tumor_subsite, therapy, sample) %>% 
  unite(sample_id_x, patient_id, sort_short_x, tumor_subsite, therapy) %>% 
  arrange(sample_id_x)

comp_tbl_sample_sort %>% 
  mutate(sort_short_x = cell_sort) %>% 
  unite(sample_id_x, patient_id, sort_short_x, tumor_subsite_rx) %>% 
  select(sample_id_x, cluster_label_sub, n, nrel, nsum) %>% 
  left_join(sample_id_x_tbl, by = "sample_id_x") %>% 
  write_tsv(paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_subtype_compositions.tsv"))
ybreaks <- c(500, 1000, 2000, 4000, 6000, 8000, 10000, 15000, 20000)

max_cells_per_sample <- max(comp_tbl_sample_sort$nsum)
ymaxn <- ybreaks[max_cells_per_sample < ybreaks][1]

comp_plot_wrapper <- function(y = "nrel", switch = NULL) {
  if (y == "nrel") ylab <- paste0("Fraction\nof cells [%]")
  if (y == "n") ylab <- paste0("Number\nof cells")
  p <- ggplot(comp_tbl_sample_sort, 
              aes_string("tumor_subsite_rx", y, fill = "cluster_label_sub")) + 
    facet_grid(~patient_id, space = "free", scales = "free", switch = switch) +
    coord_cartesian(clip = "off") + 
    scale_fill_manual(values = clrs$cluster_label_sub[[coi]]) + 
    theme(axis.text.x = element_blank(),
          axis.title.y = element_text(angle = 0, vjust = 0.5, hjust = 0.5, 
                                      margin = margin(0, -0.4, 0, 0, unit = "npc")),
          axis.ticks.x = element_blank(),
          axis.title.x = element_blank(),
          axis.line.x = element_blank(),
          strip.text.y = element_blank(),
          strip.text.x = element_blank(),
          strip.background.y = element_blank(),
          strip.background.x = element_blank(),
          plot.margin = margin(0, 0, 0, 0, "npc")) + 
    labs(x = "",
         y = ylab) +
    guides(fill = FALSE)
  if (y == "nrel") p <- p + 
    geom_bar(stat = "identity") +
    scale_y_continuous(expand = c(0, 0), 
                       breaks = c(0, 50, 100), 
                       labels = c("0", "50", "100"))
  if (y == "n") p <- p + 
    geom_bar(stat = "identity", position = position_stack()) +
    scale_y_continuous(expand = c(0, 0), 
                       breaks = c(0, ymaxn/2, ymaxn),
                       limits = c(0, ymaxn),
                       labels = c("", ymaxn/2, ymaxn)) +
    expand_limits(y = c(0, ymaxn)) +
    theme(panel.grid.major.y = element_line(linetype = 1, color = "grey90", size = 0.5))
  return(p)
} 

common_label_layers <- list(
  geom_tile(color = "white", size = 0.15),
  theme(axis.text.x = element_blank(),
        axis.ticks = element_blank(),
        axis.title.x = element_blank(),
        axis.line.x = element_blank(),
        strip.text = element_blank(),
        strip.background = element_blank(),
        plot.margin = margin(0, 0, 0, 0, "npc")),
  scale_y_discrete(expand = c(0, 0)),
  labs(y = ""),
  guides(fill = FALSE),
  facet_grid(~patient_id, 
             space = "free", scales = "free")
)

comp_label_site <- ggplot(distinct(comp_tbl_sample_sort, tumor_subsite_rx, therapy, tumor_supersite, label_supersite, patient_id), 
       aes(tumor_subsite_rx, label_supersite, 
           fill = tumor_supersite)) + 
  scale_fill_manual(values = clrs$tumor_supersite) +
  common_label_layers

comp_label_rx <- ggplot(distinct(comp_tbl_sample_sort, tumor_subsite_rx, therapy, tumor_supersite, label_therapy, consensus_signature, patient_id), 
       aes(tumor_subsite_rx, label_therapy, 
           fill = therapy)) + 
  scale_fill_manual(values = c(`post-Rx` = "gold3", `pre-Rx` = "steelblue")) +
  common_label_layers

comp_label_mutsig <- ggplot(distinct(comp_tbl_sample_sort, tumor_subsite_rx, therapy, tumor_supersite, label_mutsig, consensus_signature, patient_id), 
       aes(tumor_subsite_rx, label_mutsig, 
           fill = consensus_signature)) + 
  scale_fill_manual(values = clrs$consensus_signature) +
  common_label_layers

patient_label_tbl <- distinct(comp_tbl_sample_sort, patient_id, .keep_all = T)

comp_label_patient_id <- ggplot(patient_label_tbl, aes(tumor_subsite_rx, label_patient_id)) + 
  scale_fill_manual(values = clrs$consensus_signature) +
  geom_text(aes(tumor_subsite_rx, label_patient_id, label = patient_id)) +
  facet_grid(~patient_id, 
             space = "free", scales = "free") +
  coord_cartesian(clip = "off") + 
  theme_void() +
  theme(strip.text = element_blank(),
        strip.background = element_blank(),
        plot.margin = margin(0, 0, 0, 0, "npc"))

hist_plot_wrapper <- function(x = "nrel") {
  if (x == "nrel") {
    xlab <- paste0("Fraction of cells [%]")
    bw <- 5
  }
  if (x == "log10n") {
    xlab <- paste0("Number of cells")
    bw <- 0.2
  }
  p <- ggplot(comp_tbl_sample_sort) +
    ggridges::geom_density_ridges(
      aes_string(x, "cluster_label_sub", fill = "cluster_label_sub"), color = "black",
      stat = "binline", binwidth = bw, scale = 3) +
    facet_grid(label_supersite~., 
               space = "free", scales = "free") +
    scale_fill_manual(values = clrs$cluster_label_sub[[coi]]) + 
    theme(axis.text.y = element_blank(),
          axis.ticks.y = element_blank(),
          axis.title.y = element_blank(),
          axis.line.y = element_blank(),
          strip.text = element_blank(),
          strip.background = element_blank(),
          plot.margin = margin(0, 0, 0, 0, "npc")) +
    scale_x_continuous(expand = c(0, 0)) +
    scale_y_discrete(expand = c(0, 0)) +
    guides(fill = FALSE) +
    labs(x = xlab)
  if (x == "log10n") p <- p + expand_limits(x = c(0, 3)) + 
    scale_x_continuous(expand = c(0, 0), 
                       labels = function(x) parse(text = paste("10^", x)))
  return(p)
}

pcomp1 <- comp_plot_wrapper("n")
pcomp2 <- comp_plot_wrapper("nrel")

pcomp_grid <- plot_grid(comp_label_patient_id, 
                        pcomp1, pcomp2, 
                        comp_label_site, comp_label_rx, comp_label_mutsig,
                        ncol = 1, align = "v", 
                        rel_heights = c(0.15, 0.33, 0.33, 0.06, 0.06, 0.06))

phist1 <- hist_plot_wrapper("log10n")

pcomp_hist_grid <- ggdraw() +
  draw_plot(pcomp_grid, x = 0.01, y = 0, width = 0.85, height = 1) +
  draw_plot(phist1, x = 0.87, y = 0.05, width = 0.12, height = 0.8)

pcomp_hist_grid

# ggsave(paste0("_fig/02_composition_v6_",coi,".pdf"), pcomp_hist_grid, width = 10, height = 2)

4.3.3 site specific cluster enrichment

comp_tbl_z <- comp_tbl_sample_sort %>% 
  filter(therapy == "pre-Rx",
         !(tumor_supersite %in% c("Ascites", "Other"))) %>% 
  group_by(patient_id, cluster_label_sub) %>% 
  arrange(patient_id, cluster_label_sub, nrel) %>% 
  mutate(rank = row_number(nrel),
         z_rank = scales::rescale(rank)) %>% 
  mutate(mean_nrel = mean(nrel, na.rm = T),
         sd_nrel = sd(nrel, na.rm = T),
         z_nrel = (nrel - mean_nrel) / sd_nrel) %>% 
  ungroup()

ggplot(comp_tbl_z) +
  geom_boxplot(aes(tumor_supersite, z_nrel, color = tumor_supersite), 
               outlier.shape = NA) +
  geom_point(aes(tumor_supersite, z_nrel, color = tumor_supersite), position = position_jitter(), size = 0.1) +
  facet_grid(~cluster_label_sub, scales = "free_x", space = "free_x") +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5),
        aspect.ratio = 5) +
  scale_color_manual(values = clrs$tumor_supersite)

ggplot(comp_tbl_z) +
  geom_boxplot(aes(tumor_supersite, z_rank, color = tumor_supersite), 
               outlier.shape = NA) +
  geom_point(aes(tumor_supersite, z_rank, color = tumor_supersite), position = position_jitter(), size = 0.1) +
  facet_grid(~cluster_label_sub, scales = "free_x", space = "free_x") +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5),
        aspect.ratio = 5) +
  scale_color_manual(values = clrs$tumor_supersite)

4.4 CD8 T cell diffusion maps

## subset CD8 cells
# seu_obj_cd8 <- seu_obj_sub_sub %>%
#   subset(subset = cluster_label_sub %in% grep("^CD8.T", unique(seu_obj_sub_sub$cluster_label_sub), value = T)) %>%
#   RunUMAP(dims = 1:50, reduction = "harmony", reduction.name = "umapharmony", reduction.key = "umapharmony_")
# 
# write_rds(seu_obj_cd8, "/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/CD8.T_processed.rds")

seu_obj_cd8 <- read_rds("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/CD8.T_processed.rds")

seu_obj_cd8$cluster_label_sub <- setNames(cluster_label_extended[seu_obj_cd8$cluster_extended], seu_obj_cd8$cell_id)

write_rds(seu_obj_cd8, "/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/CD8.T_processed.rds")
# root_cell <- "SPECTRUM-OV-041_S1_CD45P_RIGHT_FALLOPIAN_TUBE_ACCCTCATCCAAGCAT"
# seu_obj_cd8_sub <- subset(seu_obj_cd8, cells = c(root_cell, colnames(seu_obj_cd8)[colnames(seu_obj_cd8)!=root_cell][-1]))
# 
# dc_obj <- DiffusionMap(seu_obj_cd8_sub@reductions$harmony@cell.embeddings, k = 100)
# dc_mat <- dc_obj@eigenvectors
# colnames(dc_mat) <- paste0("DC_", 1:ncol(dc_mat))
# seu_obj_cd8_sub[["DC"]] <- CreateDimReducObject(embeddings = dc_mat, key = "DC_", assay = DefaultAssay(seu_obj_cd8_sub))
# 
# write_rds(seu_obj_cd8_sub, "/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/CD8.T_processed_filtered.rds")
# 
# dpt_obj <- DPT(dc_obj)
# write_rds(dpt_obj, "/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/CD8.T_sub_dpt.rds")
# seu_obj_cd8_sub$DPT1 <- dpt_obj$DPT1
# write_rds(seu_obj_cd8_sub, "/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/CD8.T_processed_filtered.rds")

seu_obj_cd8_sub <- read_rds("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/CD8.T_processed_filtered.rds")

seu_obj_cd8_sub$cluster_label_sub <- setNames(cluster_label_extended[seu_obj_cd8_sub$cluster_extended], seu_obj_cd8_sub$cell_id)

write_rds(seu_obj_cd8_sub, "/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/CD8.T_processed_filtered.rds")
# plot_data_sub_sub <- as_tibble(FetchData(seu_obj_sub_sub, c(myfeatures, "cluster_label"))) %>% 
#   left_join(meta_tbl, by = "sample") %>% 
#   mutate(patient_id = str_remove_all(patient_id, "SPECTRUM-OV-"),
#          tumor_supersite = ordered(tumor_supersite, levels = names(clrs$tumor_supersite))) %>% 
#   mutate(cell_id = colnames(seu_obj_sub_sub),
#          cluster_label = ordered(cluster_label, levels = my_subtypes),
#          )
#   
# if (cell_sort == "CD45+") {
# plot_data_sub_sub <- filter(plot_data_sub_sub, sort_parameters != "singlet, live, CD45-", !is.na(tumor_supersite))
# }
# 
# if (cell_sort == "CD45-") {
# plot_data_sub_sub <- filter(plot_data_sub_sub, sort_parameters != "singlet, live, CD45+", !is.na(tumor_supersite))
# }
# 
# ggplot(plot_data_sub_sub, aes(umapharmony_1, umapharmony_2, color = cluster_label)) + 
#   common_layers_disc +
#   ggtitle("Sub cluster") +
#   #facet_wrap(~cluster_label) +
#   scale_color_manual(values = clrs$cluster_label[[coi]])
# 
# ggplot(plot_data_sub_sub, aes(umapharmony_1, umapharmony_2, color = patient_id_short)) + 
#   common_layers_disc +
#   ggtitle("Patient") +
#   # facet_wrap(~therapy) +
#   scale_color_manual(values = clrs$patient_id_short)
# 
# ggplot(plot_data_sub_sub, aes(umapharmony_1, umapharmony_2, color = tumor_supersite)) + 
#   # geom_point(aes(umapharmony_1, umapharmony_2), 
#   #            color = "grey90", size = 0.01, 
#   #            data = select(plot_data_sub_sub, -tumor_supersite)) +
#   common_layers_disc +
#   ggtitle("Site") +
#   # facet_wrap(~therapy) +
#   scale_color_manual(values = clrs$tumor_supersite)
# 
# write_tsv(select(plot_data_sub_sub, cell_id, everything(), -umapharmony_1, -umapharmony_2, -contains("RNA_")), paste0("/work/shah/uhlitzf/data/SPECTRUM/freeze/v6/", coi, "_embedding_sub.tsv"))

4.4.1 UMAP & DC

plot_data_cd8_sub <- as_tibble(FetchData(seu_obj_cd8_sub, c(myfeatures, "cluster_label_sub", "DC_1", "DC_2"))) %>% 
  left_join(meta_tbl, by = "sample") %>% 
  mutate(patient_id = str_remove_all(patient_id, "SPECTRUM-OV-"),
         tumor_supersite = ordered(tumor_supersite, levels = names(clrs$tumor_supersite))) %>% 
  mutate(cell_id = colnames(seu_obj_cd8_sub),
         cluster_label_sub = ordered(cluster_label_sub, levels = names(clrs$cluster_label_sub[[coi]])))
  
if (cell_sort == "CD45+") {
plot_data_cd8_sub <- filter(plot_data_cd8_sub, sort_parameters != "singlet, live, CD45-", !is.na(tumor_supersite))
}

if (cell_sort == "CD45-") {
plot_data_cd8_sub <- filter(plot_data_cd8_sub, sort_parameters != "singlet, live, CD45+", !is.na(tumor_supersite))
}

ggplot(plot_data_cd8_sub, aes(umapharmony_1, umapharmony_2, color = cluster_label_sub)) + 
  common_layers_disc +
  ggtitle("Sub cluster") +
  #facet_wrap(~cluster_label_sub) +
  scale_color_manual(values = clrs$cluster_label_sub[[coi]])

ggplot(plot_data_cd8_sub, aes(umapharmony_1, umapharmony_2, color = patient_id_short)) + 
  common_layers_disc +
  ggtitle("Patient") +
  # facet_wrap(~therapy) +
  scale_color_manual(values = clrs$patient_id_short)

ggplot(plot_data_cd8_sub, aes(umapharmony_1, umapharmony_2, color = tumor_supersite)) + 
  common_layers_disc +
  ggtitle("Site") +
  # facet_wrap(~therapy) +
  scale_color_manual(values = clrs$tumor_supersite)

ggplot(plot_data_cd8_sub, aes(DC_1, DC_2, color = cluster_label_sub)) + 
  common_layers_disc +
  ggtitle("Sub cluster") +
  #facet_wrap(~cluster_label_sub) +
  scale_color_manual(values = clrs$cluster_label_sub[[coi]])

ggplot(plot_data_cd8_sub, aes(DC_1, DC_2, color = patient_id_short)) + 
  common_layers_disc +
  ggtitle("Patient") +
  # facet_wrap(~therapy) +
  scale_color_manual(values = clrs$patient_id_short)

ggplot(plot_data_cd8_sub, aes(DC_1, DC_2, color = tumor_supersite)) + 
  common_layers_disc +
  ggtitle("Site") +
  # facet_wrap(~therapy) +
  scale_color_manual(values = clrs$tumor_supersite)

5 session info

devtools::session_info()
## ─ Session info ───────────────────────────────────────────────────────────────
##  setting  value                       
##  version  R version 3.6.2 (2019-12-12)
##  os       Debian GNU/Linux 10 (buster)
##  system   x86_64, linux-gnu           
##  ui       X11                         
##  language (EN)                        
##  collate  en_US.UTF-8                 
##  ctype    en_US.UTF-8                 
##  tz       Etc/UTC                     
##  date     2021-02-25                  
## 
## ─ Packages ───────────────────────────────────────────────────────────────────
##  package              * version    date       lib
##  abind                  1.4-5      2016-07-21 [2]
##  ape                    5.3        2019-03-17 [2]
##  assertthat             0.2.1      2019-03-21 [2]
##  backports              1.1.10     2020-09-15 [1]
##  bibtex                 0.4.2.2    2020-01-02 [2]
##  Biobase                2.46.0     2019-10-29 [2]
##  BiocGenerics           0.32.0     2019-10-29 [2]
##  BiocParallel           1.20.1     2019-12-21 [2]
##  bitops                 1.0-6      2013-08-17 [2]
##  boot                   1.3-24     2019-12-20 [3]
##  broom                  0.7.2      2020-10-20 [1]
##  callr                  3.4.2      2020-02-12 [1]
##  car                    3.0-8      2020-05-21 [1]
##  carData                3.0-4      2020-05-22 [1]
##  caTools                1.17.1.4   2020-01-13 [2]
##  cellranger             1.1.0      2016-07-27 [2]
##  class                  7.3-15     2019-01-01 [3]
##  cli                    2.0.2      2020-02-28 [1]
##  cluster                2.1.0      2019-06-19 [3]
##  codetools              0.2-16     2018-12-24 [3]
##  colorblindr          * 0.1.0      2020-01-13 [2]
##  colorspace           * 1.4-2      2019-12-29 [2]
##  cowplot              * 1.0.0      2019-07-11 [2]
##  crayon                 1.3.4      2017-09-16 [1]
##  curl                   4.3        2019-12-02 [2]
##  data.table             1.12.8     2019-12-09 [2]
##  DBI                    1.1.0      2019-12-15 [2]
##  dbplyr                 2.0.0      2020-11-03 [1]
##  DelayedArray           0.12.2     2020-01-06 [2]
##  DEoptimR               1.0-8      2016-11-19 [1]
##  desc                   1.2.0      2018-05-01 [2]
##  destiny              * 3.0.1      2020-01-16 [1]
##  devtools               2.2.1      2019-09-24 [2]
##  digest                 0.6.25     2020-02-23 [1]
##  dplyr                * 1.0.2      2020-08-18 [1]
##  e1071                  1.7-3      2019-11-26 [1]
##  ellipsis               0.3.1      2020-05-15 [1]
##  evaluate               0.14       2019-05-28 [2]
##  fansi                  0.4.1      2020-01-08 [2]
##  farver                 2.0.3      2020-01-16 [1]
##  fitdistrplus           1.0-14     2019-01-23 [2]
##  forcats              * 0.5.0      2020-03-01 [1]
##  foreign                0.8-74     2019-12-26 [3]
##  formattable            0.2.0.1    2016-08-05 [1]
##  fs                     1.5.0      2020-07-31 [1]
##  future                 1.15.1     2019-11-25 [2]
##  future.apply           1.4.0      2020-01-07 [2]
##  gbRd                   0.4-11     2012-10-01 [2]
##  gdata                  2.18.0     2017-06-06 [2]
##  generics               0.0.2      2018-11-29 [2]
##  GenomeInfoDb           1.22.0     2019-10-29 [2]
##  GenomeInfoDbData       1.2.2      2020-01-14 [2]
##  GenomicRanges          1.38.0     2019-10-29 [2]
##  ggplot.multistats      1.0.0      2019-10-28 [1]
##  ggplot2              * 3.3.2      2020-06-19 [1]
##  ggrepel                0.8.1      2019-05-07 [2]
##  ggridges               0.5.2      2020-01-12 [2]
##  ggthemes               4.2.0      2019-05-13 [1]
##  globals                0.12.5     2019-12-07 [2]
##  glue                   1.3.2      2020-03-12 [1]
##  gplots                 3.0.1.2    2020-01-11 [2]
##  gridExtra              2.3        2017-09-09 [2]
##  gtable                 0.3.0      2019-03-25 [2]
##  gtools                 3.8.1      2018-06-26 [2]
##  haven                  2.3.1      2020-06-01 [1]
##  hexbin                 1.28.0     2019-11-11 [2]
##  hms                    0.5.3      2020-01-08 [1]
##  htmltools              0.5.1.1    2021-01-22 [1]
##  htmlwidgets            1.5.1      2019-10-08 [2]
##  httr                   1.4.2      2020-07-20 [1]
##  ica                    1.0-2      2018-05-24 [2]
##  igraph                 1.2.5      2020-03-19 [1]
##  IRanges                2.20.2     2020-01-13 [2]
##  irlba                  2.3.3      2019-02-05 [2]
##  jsonlite               1.7.1      2020-09-07 [1]
##  KernSmooth             2.23-16    2019-10-15 [3]
##  knitr                  1.26       2019-11-12 [2]
##  labeling               0.3        2014-08-23 [2]
##  laeken                 0.5.1      2020-02-05 [1]
##  lattice                0.20-38    2018-11-04 [3]
##  lazyeval               0.2.2      2019-03-15 [2]
##  leiden                 0.3.1      2019-07-23 [2]
##  lifecycle              0.2.0      2020-03-06 [1]
##  listenv                0.8.0      2019-12-05 [2]
##  lmtest                 0.9-37     2019-04-30 [2]
##  lsei                   1.2-0      2017-10-23 [2]
##  lubridate              1.7.9.2    2020-11-13 [1]
##  magrittr             * 2.0.1      2020-11-17 [1]
##  MASS                   7.3-51.5   2019-12-20 [3]
##  Matrix                 1.2-18     2019-11-27 [3]
##  matrixStats            0.56.0     2020-03-13 [1]
##  memoise                1.1.0      2017-04-21 [2]
##  metap                  1.2        2019-12-08 [2]
##  mnormt                 1.5-5      2016-10-15 [2]
##  modelr                 0.1.8      2020-05-19 [1]
##  multcomp               1.4-12     2020-01-10 [2]
##  multtest               2.42.0     2019-10-29 [2]
##  munsell                0.5.0      2018-06-12 [2]
##  mutoss                 0.1-12     2017-12-04 [2]
##  mvtnorm                1.0-12     2020-01-09 [2]
##  nlme                   3.1-143    2019-12-10 [3]
##  nnet                   7.3-12     2016-02-02 [3]
##  npsurv                 0.4-0      2017-10-14 [2]
##  numDeriv               2016.8-1.1 2019-06-06 [2]
##  openxlsx               4.1.5      2020-05-06 [1]
##  pbapply                1.4-2      2019-08-31 [2]
##  pcaMethods             1.78.0     2019-10-29 [2]
##  pillar                 1.4.6      2020-07-10 [1]
##  pkgbuild               1.0.6      2019-10-09 [2]
##  pkgconfig              2.0.3      2019-09-22 [1]
##  pkgload                1.0.2      2018-10-29 [2]
##  plotly                 4.9.1      2019-11-07 [2]
##  plotrix                3.7-7      2019-12-05 [2]
##  plyr                   1.8.5      2019-12-10 [2]
##  png                    0.1-7      2013-12-03 [2]
##  prettyunits            1.1.1      2020-01-24 [1]
##  processx               3.4.2      2020-02-09 [1]
##  progeny              * 1.11.3     2020-10-22 [1]
##  proxy                  0.4-24     2020-04-25 [1]
##  ps                     1.3.2      2020-02-13 [1]
##  purrr                * 0.3.4      2020-04-17 [1]
##  R.methodsS3            1.7.1      2016-02-16 [2]
##  R.oo                   1.23.0     2019-11-03 [2]
##  R.utils                2.9.2      2019-12-08 [2]
##  R6                     2.4.1      2019-11-12 [1]
##  ranger                 0.12.1     2020-01-10 [1]
##  RANN                   2.6.1      2019-01-08 [2]
##  rappdirs               0.3.1      2016-03-28 [2]
##  RColorBrewer           1.1-2      2014-12-07 [2]
##  Rcpp                   1.0.4      2020-03-17 [1]
##  RcppAnnoy              0.0.16     2020-03-08 [1]
##  RcppEigen              0.3.3.7.0  2019-11-16 [2]
##  RcppHNSW               0.2.0      2019-09-20 [2]
##  RcppParallel           4.4.4      2019-09-27 [2]
##  RCurl                  1.98-1.1   2020-01-19 [1]
##  Rdpack                 0.11-1     2019-12-14 [2]
##  readr                * 1.4.0      2020-10-05 [1]
##  readxl               * 1.3.1      2019-03-13 [2]
##  rematch                1.0.1      2016-04-21 [2]
##  remotes                2.1.0      2019-06-24 [2]
##  reprex                 0.3.0      2019-05-16 [2]
##  reshape2               1.4.3      2017-12-11 [2]
##  reticulate             1.14       2019-12-17 [2]
##  rio                    0.5.16     2018-11-26 [1]
##  rlang                  0.4.8      2020-10-08 [1]
##  rmarkdown              2.0        2019-12-12 [2]
##  robustbase             0.93-6     2020-03-23 [1]
##  ROCR                   1.0-7      2015-03-26 [2]
##  rprojroot              1.3-2      2018-01-03 [2]
##  RSpectra               0.16-0     2019-12-01 [2]
##  rstudioapi             0.11       2020-02-07 [1]
##  rsvd                   1.0.3      2020-02-17 [1]
##  Rtsne                  0.15       2018-11-10 [2]
##  rvest                  0.3.6      2020-07-25 [1]
##  S4Vectors              0.24.2     2020-01-13 [2]
##  sandwich               2.5-1      2019-04-06 [2]
##  scales                 1.1.0      2019-11-18 [2]
##  scatterplot3d          0.3-41     2018-03-14 [1]
##  sctransform            0.2.1      2019-12-17 [2]
##  SDMTools               1.1-221.2  2019-11-30 [2]
##  sessioninfo            1.1.1      2018-11-05 [2]
##  Seurat               * 3.1.2      2019-12-12 [2]
##  SingleCellExperiment   1.8.0      2019-10-29 [2]
##  smoother               1.1        2015-04-16 [1]
##  sn                     1.5-4      2019-05-14 [2]
##  sp                     1.4-2      2020-05-20 [1]
##  stringi                1.5.3      2020-09-09 [1]
##  stringr              * 1.4.0      2019-02-10 [1]
##  SummarizedExperiment   1.16.1     2019-12-19 [2]
##  survival               3.1-8      2019-12-03 [3]
##  testthat               2.3.2      2020-03-02 [1]
##  TFisher                0.2.0      2018-03-21 [2]
##  TH.data                1.0-10     2019-01-21 [2]
##  tibble               * 3.0.4      2020-10-12 [1]
##  tidyr                * 1.1.2      2020-08-27 [1]
##  tidyselect             1.1.0      2020-05-11 [1]
##  tidyverse            * 1.3.0      2019-11-21 [2]
##  tsne                   0.1-3      2016-07-15 [2]
##  TTR                    0.23-6     2019-12-15 [1]
##  usethis                1.5.1      2019-07-04 [2]
##  uwot                   0.1.5      2019-12-04 [2]
##  vcd                    1.4-7      2020-04-02 [1]
##  vctrs                  0.3.5      2020-11-17 [1]
##  VIM                    6.0.0      2020-05-08 [1]
##  viridis              * 0.5.1      2018-03-29 [2]
##  viridisLite          * 0.3.0      2018-02-01 [2]
##  withr                  2.3.0      2020-09-22 [1]
##  xfun                   0.12       2020-01-13 [2]
##  xml2                   1.3.2      2020-04-23 [1]
##  xts                    0.12-0     2020-01-19 [1]
##  XVector                0.26.0     2019-10-29 [2]
##  yaml                   2.2.1      2020-02-01 [1]
##  zip                    2.0.4      2019-09-01 [1]
##  zlibbioc               1.32.0     2019-10-29 [2]
##  zoo                    1.8-7      2020-01-10 [2]
##  source                                 
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  Bioconductor                           
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Github (clauswilke/colorblindr@1ac3d4d)
##  R-Forge (R 3.6.2)                      
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  Bioconductor                           
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Github (saezlab/progeny@94bea1c)       
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.3)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.3)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
##  CRAN (R 3.6.2)                         
##  Bioconductor                           
##  CRAN (R 3.6.2)                         
## 
## [1] /home/uhlitzf/R/lib
## [2] /usr/local/lib/R/site-library
## [3] /usr/local/lib/R/library